EIGENFUNCTIONS AND VERY SINGULAR SIMILARITY 
SOLUTIONS OF ODD-ORDER NONLINEAR DISPERSION PDES 



R.S. FERNANDES AND V.A. GALAKTIONOV 

Abstract. Asymptotic properties of solutions of the nonlinear dispersion equations 

(0.1) u t = {\u\ n u) xxx and u t = {\u\ n u) xxx - \u\ p ~ 1 u in Ixl + , 

where n > and p > n + 1 are fixed exponents, and their higher-odd-order analogies 
are studied. The global in time similarity solutions, which lead to "nonlinear eigenfunc- 
tions" of the rescaled ODEs, are constructed. The basic mathematical tools include a 
"homotopy-deformation" approach, where the limit n —> + in the first equation in (|0.1[) 
turns out to be fruitful by reducing the problems to the linear one at n — 0, i.e., 

v XXXJ 

for which Hermitian spectral theory was developed earlier in [4] , as well as other nonlinear 
operator and numerical methods. The nonlinear limit n — > +oo, with the limit PDE 

(signw) t = v xxx , in terms of the variable v = \u\ n u, 

admitting three almost explicit nonlinear eigenfunctions, has also been described. 

For the second equation in (jO.ip . very singular similarity solutions (VSSs) are con- 
structed. A "nonlinear bifurcation" phenomenon at critical values {p — pi(n)}i>o of the 
absorption exponents is discussed. In fact, regardless their wide spread in applications 
and clear significance, such third-order nonlinear dispersion PDEs are quite poorly repre- 
sented in general PDE theory, unlike their direct "neighbouring" second and fourth-order 
quasilinear parabolic equations (the so-called porous medium equations, PMEs) 

u t = (\u\ n u) xx (the PME-2) and u t = -{\u\ n u) xxxx (the PME-4). 

The main goal of the present paper is to essentially fill such a gap. 



1. INTRODUCTION: NONLINEAR DISPERSION PDEs 

1.1. NDEs: global smooth similarity patterns. In the present paper, we study as- 
ymptotic properties of nonlinear dispersion equations (NDEs) of the following form: 

(1.1) u t = (-l) k+l Df + \\u\ n u)-\u\ p - 1 u in Kxl + , k = l,2,..., 
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where n > and p > n + 1 are fixed exponents. This continues the study began in [1] 
for n — 0, i.e., for the semilinear dispersion equation. A continuous "homotopy" path 
n — > + connecting the present study with that in pE] turns out to be rather effective. 

Let us clarify from the beginning a general role and the main purpose of the present 
study in the theory of higher-order NDEs such as (II. 1ft : 

(i) We will study here a countable family of sufficiently smooth continuous similarity 
solutions of VSS (very singular) type, which are defined for all t > 0. 

(ii) In addition, (II. ip . as an equation with nonlinear dispersion mechanism, contains 
and describes other key singularity phenomena such as a complicated formation of various 
shock and smoother rarefaction waves, which appear from discontinuous data, as well as a 
general nonuniqueness of "entropy solutions" after shocks. We do not touch these difficult, 
even mathematically obscure, phenomena and refer to [SJ [T31 [5] for further details. 

Overall, it may be said that the smooth similarity behaviour and asymptotic patterns 
studied here occur in the NDE (II. ip for large times, when the shock wave influence has 
already been settled down via evolution, and becomes negligible. 

Concerning local existence and uniqueness theory, including "smoothing results" for 
NDEs (i.e., for solutions without shocks and stronger singularities), see references and 
results in [8]. Since we are dealing with some special exact similarity solutions of (II. ip . 
we do not use any advanced results of local and/or global regularity and any shock-entropy 
theory here, and always tackle continuous solutions. 

Various applications of nonlinear dispersion equations are characterized in detail in JT3] , 
so we will minimize extra references to previous physical and more formal investigations 
of such PDE models, and concentrate on their quite unusual mathematical aspects of 
our interest. It is worth mentioning here a couple of remarkable examples of NDEs from 
integrable PDE theory. The first one is the third-order Harry Dym equation 

(1.2) u t = u 3 u xxx , 

a most exotic integrable soliton equation; see [HI § 4.7] for survey and references therein. 
The second one is the fifth-order Kawamoto equation [T6] 

U U X xxxx ~\~ 5 U U x U X xxx 10 11 U xx U X xxi 

which has higher-degree algebraic terms. Shock waves for the pure fifth-order NDEs 

iH ^ Hxxxxx and itj (|^| ^) xxxxxi 

are studied in [H] . Let us discuss further necessary aspects of NDEs and their standing in 
general PDE theory. 

1.2. Compactons in NDEs: compactly supported travelling waves. Consider the 
higher odd-order nonlinear dispersion equation with another divergent lower-order term: 

(1.3) u t = {-l) k+1 Dl k+ \\u\ n u) + (|t*|"*) x . 

Equation ( II. 3p is a generalization of the third-order Rosenau-Hyman (RH) equation 
(1-4) u t = {u 2 ) xxx + (u 2 ) x , 
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which models the effect of nonlinear dispersion in the pattern formation of liquid drops 
(see [2T]) . For n = k = 1, (ll.3p is the RH equation in classes of nonnegative solutions. 

It is well known that the RH equation ( j 1.4ft possesses explicit moving compactly sup- 
ported soliton-type solutions known as compactons. Compactons have the same structure 
as travelling wave solutions, given by 

u c (x,t) — f(z), where z = x — At. 

So looking at compacton solutions for (II .3p . on substitution we have that 

-xr = (-i) k+i D? + \\frf) + a/ivy. 

After integrating once with zero constant (i.e., zero "flux"), we find that f(z) satisfies 

(1.5) - X f = (-l) k+1 D 2 z k (\f\ n f) + \f\ n f. 

For k = 1, the ODE (11.51) can be solved explicitly (see the expression (11.91) below), while 
for k > 2, this is a difficult variational problem, which admits various countable families 
of compactly supported oscillatory solutions f(z) of changing sign, [12J. 

Whilst we note that compacton solutions may be found for nonlinear dispersion equa- 
tions, for non-conservative and non-divergent NDEs such as (11.11) . TW solutions may be 
irrelevant for classes of bounded solutions of the Cauchy problem. Therefore, instead we 
look to find other more complicated similarity solutions. 

1.3. A relation to blow-up in reaction-diffusion theory. Surprisingly, the NDE 
(II. 3p is somehow related to the parabolic even-order equations of reaction-diffusion type: 

(1.6) u t = {-l) k+1 D 2 x k {\u\ n u) + \u\ n u (u t = {u n+1 ) xx + u n+1 for k = 1, u > 0). 
Equation (II. 6p admits blow-up self-similar solutions of the separate form 

(1.7) u(x,t) = (T-t)-$f(x), 

where T is the finite blow-up time. This self-similarity reduces the PDE to the ODE for 
the similarity profile / easily obtained by substituting (I1.7P to f II . 3 p : 

(1.8) (-i)* +l £>?(i/r/) + i/r/=i/- 

For k = 1, equation (ll.Sp possesses the explicit compactly supported solution 



:i-9) m 



" 2(n+1) cos 2 ( nx VI" for \x\ < ^ (n+1) 

_n(n+2) LUb \2{n+l))\ iUI n ■ 

for \x\ > ^±^. 



Thus, (II. 7p . (II. 9p yields the so-called standing wave blow-up solution (S-regime of blow- 
up) of (II. 3p . which have compact support for all times t 6 (0,T). This exact Zmitrenko- 
Kurdyumov blow-up solution has been known since the middle of 1970s; see details in [22| 
Ch. 4]. As we have mentioned, for k > 2, (I1.8P cannot be solved explicitly, but the ODE is 
shown to admit a countable set of compactly supported solutions obtained by variational 
Lusternik-Schnirel'man and fibering methods, [T2] . 
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Thus, comparing the ODEs (I1.8P and (|1.5|) . we see that these coincide provided the 
compacton velocity A is given by 

(1.10) A = -i. 

In other words, this means that some principles of compacton propagation for such NDEs 
(11 . 3 p are directly related to blow-up formation mechanisms for reaction-diffusion equa- 
tions fll.6p . Moreover, in both equations, such asymptotic structures are expected to be 
structurally stable (for k — 1 in (11.61) . this has been proved, (221 P- 260]). This is in- 
deed surprising, since both classes of nonlinear PDEs seem to be responsible for entirely 
different physical phenomena (to say nothing of their different mathematical essence). 

The lower order case of equation (11. 6p for k — 1, which is just a standard reaction- 
diffusion PDE, is fairly well understood. However, the third-order nonlinear dispersion 
equation in (II .3p . also for the minimal value k — 1, has not been studied extensively and 
some of its basic governing mathematical principles are still relatively unknown. 

1.4. Two main model NDEs and layout of the paper. Thus, we will construct 
global similarity solutions of the following two NDEs. The first one is the pure NDE 

(1.11) u t = {-l) k+1 D 2 x k + \\u\ n u) in RxR + , n> 0, 

which is studied in Sections [2HH Our goal therein is to show that (II. lip admits an 
infinite countable family of self-similar solutions governed by "nonlinear eigenfunctions" 
of a rescaled operator. Moreover, in Section I3.6[ we show that this countable family as 
n — > + can be described by eigenfunctions from Hermitian spectral theory of the non-self 
adjoint operator [U § 4] 

(1-12) B = (-ir^ + J i ri yD y + J i ri I. 

The above operator occurs after similar scaling of the linear dispersion equation (LDE) 

(1.13) u t = {-l) k+l D 2k+l u in Rxl + ; k > 1, 

so we suggest to describe nonlinear eigenfunctions of the NDE (11.111) by using linear 
theory. 

In Section HI we perform the opposite "nonlinear" limit n — > +oo in (11.111) . Namely, we 
observe that the natural change of the independent variable leads to the following PDE: 

(1.14) \u\ n u = v =^ (\v\-^v) t = {-l) k+1 Df +1 v. 
The formal limit n — > +oo leads to the following limit NDE: 

(1.15) (signv) t = (-l) k+1 D 2 x k+1 v, 

which admits analogous similarity patterns. It turns out that, at least for k — 1, where 
(11.151) takes a particularly simple form 

(1-16) (signv) t = v xxx , 

first three occurring ODEs for (11.161) admit reducing to an algebraic system, and we 
develop a geometric-algebraic approach to constructing first nonlinear eigenfunctions. 
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In Section^ we return to the study of VSSs of the NDE with absorption (11 .11) . where our 
main goal is to justify existence of the so-called p-bifurcation branches of VSS solutions, 
which appear at some critical exponents p = pi(n) > n + 1, / = 0, 1, 2, ... . 



2. Similarity solutions of the NDE ("nonlinear eigenfunctions") 

Here, we consider the pure NDE (11.111) . which is connected to (II. 3p . but now we do not 
have the convection-like term (|w| n w) x , which is negligible in our asymptotics. This non- 
linear dispersion equation may be compared with the even-order model, which represents 
the general higher-order porous medium equation (the PME-2m) 

(2.1) u t = (-l) m+1 A m (HV) in R N xi + , m > 2. 

The classic PME-2, for m — 1, appears in a number of physical applications, such as 
fluid and gas flows, heat transfer, or (nonlinear) diffusion. Other applications have been 
proposed in mathematical biology, lubrication and boundary layer theory, and various 
other fields. Papers exploring delicate asymptotics for the PME include [TU] , where extra 
references are available. For nonlinear eigenfunctions of (12.11) with m — 2, N — 1, see 
[7] and key references therein. Overall, the higher-order case m > 2 in (12.11) has been 
studied much less in the mathematical literature and represents a number of difficult open 
problems. 

2.1. Self-similar solutions: towards a "nonlinear eigenvalue problem" . Our NDE 
(11.111) admits standard similarity solutions given by 

(2.2) u gl (x,t) =t~ a f(y), where y = x/t p , 

for some unknown real parameters a and (3. After substitution (I2.2p into the NDE, we 
obtain the ODE for the rescaled profile /, 

(2.3) - at' - 1 / - f3r a - l fy = (_l)*+l t -«(n+l)-^(2fc+l) £) 2fc+l^|n^ - 

By equating powers of t, the parameter /3 can be found in terms of a and is given by 

(2.4) £ = ^TT>0 for a<± 

so that now a G M remains the only unknown nonlinear eigenvalue. 
Thus, the ODE Q takes the form 

(2.5) A(f,a) = (-l) k+1 Dl k+1 (\f\ n f) + ^f>y + af = in R. 

In order to formulate a suitable nonlinear eigenvalue problem for the pairs {ai, fi}i>o for 
the ODE (12. 5p . one should specify the "boundary" conditions at y = ±oo, which is one of 
the main goals of this paper. It turns out that, loosely speaking, such conditions can be 
formulated as follows: for some special values of the real eigenvalues {ai = a>i{n) > 0};>o, 
the corresponding nonlinear eigenfunctions 

■ fi(y) have: (i) finite left-hand interface (fi(y) = for y <^ —1), 
and (ii) "minimal oscillatory" behaviour as y — > +oo. 



Both such conditions will get proper explanations later on. Of course, at least for suf- 
ficiently large I > 2k + 1, the problem (I2.5p . (I2.6P exhibits all typical features of self- 
similarity of the second kind (the first kind corresponds to standard dimensional analysis 
of PDEs), a notion introduced by Ya.B. Zel'dovich in 1956 [21]. Unlike previous and 
known examples, in general, the eigenvalues oti{n) (and hence /3/(n)) cannot be found 
explicitly from any dimensional analysis. Thus, admitted (real) values of a = ai(n) > 
at this stage are still unknown and play a role of "nonlinear eigenvalues" . 

In particular, our goal is to show by a combination of analytic, formal, and numerical 
tools that, for any n > 0, the eigenvalue problem (12. 5p . (12. 6 p 

(2.7) has a countable set of pairs {«/, fi(y)}i>o- A(/;, ai) = 0. 

For the linear case with n = 0, the corresponding linear non-self-adjoint Hermitian 
spectral problem was developed in [H § 4]; we present and use these results below. 

2.2. Towards blow-up patterns. The present analysis admits a natural duality: since 
the NDE (ll.lip is invariant under reflections, 

the global patterns (12. 2p lead to the corresponding reflected blow-up ones: 

(2.9) u vl (x,t) = (T-t)- a f(y), y = - w ^ F , where = > . 

Therefore, the nonlinear eigenvalue problem (12. 5p . (12. 6p is assumed to describe both count- 
able families of global and blow-up patterns for the NDE (II. lip . In what follows, for 
simplicity, we will mainly use a global treatment of such asymptotic patterns. 

2.3. A homotopy path n — > + to Hermitian spectral theory: a route to count- 
ability of nonlinear eigenfunctions. For n = 0, a gives the eigenvalues {A/} in linear 
Hermitian spectral theory, [U § 4]. Indeed, for n = 0, (I2.5P reads 

(2.10) {-lf^Df + 1 / + ^ f'y + af = Bf + (a - = 0. 
Therefore, this gives the eigenvalue equation for the linear operator B in (II. 12ft : 

(2.11) B?/> = Xip, where A = —a + ^i+i- 

By Hermitian spectral theory [4], this defines a countable set of eigenvalues for n = 0: 

(2.12) ai {Q)=i±L i Z = ,l,2,.... 

The corresponding eigenfunctions are then the normalized derivatives of the rescaled 
kernel F(y) of the fundamental solution of the dispersion operator D t + (— l) k D^ +1 : 

(2-13) ^( y ) = ^23jF(y) l l>0. 

Moreover, the "adjoint" (not in the standard L 2 -sense) operator 
(2.14) B^i-l^D^-^yDy 
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has the same spectrum and eigenfunctions {ipf}, which are generalized Hermite polyno- 
mials given by (see jU § 5.2]) 

(2.15) = ^ y'+(-i) fc+1 E jrf + V 

Finally, the operator pair {B,B*} has a bi-orthonormal sets of eigenfunctions. 

As a key conclusion, we expect that there exists a uniform and pointwise convergence 
as n — > + of the nonlinear eigenfunctions to the linear ones for B defined by (I2.13p . 
We postpone explanations of main aspects of such a branching until Section 13.61 when 
we will have gained enough understanding of basic nonlinear eigenfunction theory in- 
volved. However, a rigorous proof of such a branching of eigenfunctions at n = is a hard 
open problem. Nevertheless, this n-branching approach still remains a unique analytically 
convincing fact towards existence of an infinite countable discrete set of nonlinear eigen- 
functions of the problem (12.51) . In this case, branching theory as n — > + can be developed 
along the same lines as for the PME-4 (12 .ip with m — 2 [7J § 6]. Such a construction 
then uses a specific Hermitian spectral theory created in [I] ; see Section 13 .61 

2.4. Conservation laws: explicit values of first nonlinear eigenvalues. As a first 
pleasant surprise, it turns out that some first eigenvalues for any n > can be calcu- 
lated explicitly by conservation laws for the NDE. Assuming that the solution u(x,t) is 
integrable in x over K (this is a formal assumption that can be got rid of), we have that 
(11.111) is conservative in mass, and so 

(2.16) f t Ju(x,t)dx = 0. 

R 

For similarity solutions (12.21) . we have that 

Ju(x,t)dx = t^ a Jf(y) dy. 



1 1 

l 



This satisfies (12.161) . provided that, in addition to (12.41) . 
(2.17) (3-a = =^ ao(n) - (2fc+1)+n , 

for non-zero rescaled mass J f =fi 0. So, on substitution into (12. 51) . 

(-ir^d/IV) + fy + / = 0. 

Integrating once with the zero constant (a zero-flux condition), we end up with the ODE 
(2-18) (-1)^^(1/1-/) + fy = 0. 

Note that, for n = 0, we have exactly the linear ODE (see [1]) 

(2.19) + Fy = for ye R. 
For convenience, we use in (12. 18f) the natural substitution 

(2.20) Y=\f\ n f =^ f=\Y\-^Y, 



in order to remove nonlinearities in the highest differential. Substitution yields 
(2-21) (—l) k+1 Dy k Y + y\Y\-^Y = 0. 

Similarly, we have conservation of the first moment, with 

Jxu(x,t)dx = t 2 ^ a Jyf(y) dy. 

M R 

Hence, we have that 

(2.22) 2/3 -a = =► ai (n) = 

This then gives the ODE 

(2-23) (-i)^f +HI/IV) + prk ^ + ptW / = o- 

However, we cannot simply integrate this equation, as we could before, to reduce the 
order of the ODE. Instead, we multiply f l2.23|) by y, so that 

(-i) k+1 D?+\\frf)y + j^- n f'y 2 + p+f^ fy = o, 

and now it is possible to integrate by parts, to obtain 
(2.24) (-l) k+1 Df(\m)y + (-l) k Df- l (\f\ n f) + fy 2 = 0. 

We also look at conservation of the second moment, 

fx 2 u(x, t) dx = t 3 P- a jy 2 f{y) dy. 

R R 

This gives 3/3 — a = 0, so that, by (12.41) . 

(2-25) ^(nHjg^ and (-l) k+1 D^(\f\ »f) + /*„ + ^3 — / = 0. 

Similarly, we multiply by and integrate to reduce the order, to obtain the ODE 
{ _ 1)k+ i D 2 k{lfl n f)y 2 + 2 {-lf Df -\\f\ n f)y 

(2 ' 26) +2(-i) fc+ ^f - 2 (|/|V) + (^fe fy' = o. 

These three conservation laws in particular are important, as we can explicitly find the 
first three (second-order) equations for the case k = 1 (corresponding to the first three 
nonlinear eigenvalues), but not for others. The case k — 1 is essential, since it is much 
easier to develop theory for the lower-order CcLSC, clS well as it is easier to solve numerically 
(see Section 1231) . 

In general, for arbitrary k > 1, for all I < 2k + 1, where I is the eigenvalue index as 
before, we have our moments conservation given by 

fx l u(x, t) dx = t {l+1) ?- a Jy l f(y) dy. 



Therefore, our nonlinear eigenvalues are represented by 

(2.27) (Z + l)/3-a = =>► t* t (n) = (2fc+1 ^ +1)n for any < I < 2k + 1. 
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The corresponding ODEs, which can be integrated once for all such eigenvalues, are 

(2.28) (-l) k+l D^\\f\«f) + (2fc+1) i (t+1)w f'y + (2fc+i ;+ 1 (;+1)n / = 0. 

All of them admit reducing the order by multiplying by y l and integrating by parts I 
times. The function Y = \f\ n f then solved the semilinear PDE 

(2 29) f_i ifc+iy(2fc+i) j i v(\Y\~^YY H ^ \Y\~^Y = 

\ / ^ (2fc+l) + («+l)n y^l 1 I 1 > ^ (2fc+l)+(J+l)n 1 1 

Thus, for such NDEs, one can always obtain explicitly 2k + 1 n-branches of nonlinear 
eigenvalues (12.271) . though of course the solvability of the corresponding eigenvalue equa- 
tions remains a difficult open problem, especially for large I. However, (I2.27P shows one 
important feature of this nonlinear eigenvalue theory: the n-branches (I2.27P are global in 
n > 0, i.e., exist for all n > 0. Existence of such a global continuation, with no turning 
points, is one of the most difficult question in general nonlinear operator theory; see e.g., 
[T8| [23] . Recall that the eigenvalue problem (12. 5p is by no means variational, neither 
contains any monotone or coercive operators. 

Next, we must admit that, for any larger I > 2k+l, we cannot find ai(n) explicitly using 
conservation laws, and then, searching for n-branches of these nonlinear eigenfunctions, 
we either should rely on the above local homotopy approach as n — > + , or use advanced 
numerical methods. 

2.5. Numerical construction of nonlinear eigenfunctions. We look to use a shoot- 
ing method in order to find reliable profiles for the NDE (11.111) . First, considering solutions 
satisfying conservation of mass, we have our first "nonlinear eigenvalue" (for 1 = 0), where 
n = corresponds to the linear kernel, — F(y), i.e., A = in (12. lip . Here the rescaled 
equation is given by (12.211) . 

Consider first the simpler case k = 1, when we have the second-order equation: 
(2.30) Y" = -^\Y\-^Yy. 

Then, it can be shown that the left-hand interface is not oscillatory (in fact, the same is 
true for n = 0, i.e., for the Airy function). Hence, for small solutions of Y(y) as y — > y$, 
with some interface point yo < 0, we can approximate it by 

(2-31) Y{y) = C (y-y )%(l + o(l)) 1 

for some constant Cq > and exponent a > 0. Here we have used the standard notation 

(•)+ = max{0, •}. 

Substituting (12.311) into the ODE (I2.30p . we obtain the leading equality 

(2.32) a(a - l)C (y - y f " 2 = C^(y - y )&y , 

where we use (y — yo) to mean (y — y ) + (l + o(l)), as defined before. Hence, we must 
have from f l2T32"|) that 

2(n+l) 

a = — -, 
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with 2 < a < 4 for n > 1. From this, the constant Cq is given by 



n+l 



Co = ( 2 (n+i)(n+2)(n+3) ) " for an y interface y < 0. 
Our solution then has the expansion 

(2-33) Y(y) = L n+1 ^ in+3) )^ (V ~ Vo) 2J ^(l + o(l)) 



>2(ra+l)(n+2)(n+3) 

with the derivative expansion obtained similarly. 



as y ->• Vo 



(2-34) Y'(y) = ^( 2{n+1 ^li ){n+3) )^ (y ~ y )^(l + o(l)). 

We use the MatLab IVP solver ode 15s to plot our profiles. Taking an arbitrary initial 
(interface) point y = yo < 0, our solution and first derivative will be zero here. Since both 
initial values are zero at y = yo, we will often find the solution Y = 0, whilst trying to 
solve numerically. In order to overcome this problem, we must look at some point yo + 5 
(for small 5 ~ 10~ 3 ), close to this point, and after finding the derivative there, this is used 
as the initial condition. 

Obviously due to the nature of the expansion ( I2.34p . we must take a relatively large 
initial point (in our case we take yo = —10), in order for us to have an initial condition that 
is not negligible. Hence, the solution is a large rescaling of the "nonlinear fundamental 
solution" with the unit mass. Below are a few profiles that have been found, in which 
we have taken 5 = 10 -3 . For n ~ 0.5, the derivative Y' is very small and larger negative 
interface points must be used to find reliable profiles. This makes comparison, between 
different values of n, more difficult. 

In Figures [TH31 we show the first nonlinear eigenfunction Yo(y), with k — 1, for various 
values of n = 3,2,1, and 0.7. Note that this is precisely the profile that, as n — > + , 
must converge to the rescaled kernel F(y) = Ai(y), representing Airy's classic function 
satisfying (cf. f)2.30p for n — 0) 

(2.35) F" + | Fy = 0, J F = 1, where F E L*(R), p{y) = jj^ ^ * < £ 

where a > is a small enough constant. 

In general, for all values of I, we have, for the lower-order case of k = 1, that the same 
expansion (I2.33j) and (I2.34p hold near finite interfaces. Here I < 2k + 1 and hence, for 
k = 1, we can only take / = 0,1,2 for explicitly given nonlinear eigenvalues. 

For / = 1, we have from (I2.24p . with Y(y) = \f\ n f and k — 1, that the eigenfunction 
Yi(y) solves the ODE 

(2-36) Y »=lY-^\Y\-^Yy. 
For I = 2 and k = 1, we have from ( 12.26P that Y 2 (y) solves 
(2-37) Y » = *Y>-$Y-^\Y\-**Yy. 

We use these to plot the profiles of our NDE with / = 1 and 1 = 2, which, for n = 0, 
coincide with the derivatives F'(y) and F"{y) respectively, of the linear kernel F(y). 
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FIGURE 1. Rescaled solution Y(y) of the ODE (EHUl) for k = 1, I = 0, with n = 3. 
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FIGURE 2. Rescaled solution Y(y) of the ODE fl!T30|) for A; = 1, I = 0, with n = 2. 

Figures [5] and [6] show the "dipole-like" profiles for / = 1 as solutions of (I2.36P , while 
Figure [7] represents even more oscillatory third eigenfunction /^(y) for / = 2. Recall 
that, by a homotopy path as n — > + , this function is expected to converge to the highly 
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FIGURE 3. Rescaled solution Y(y) of the ODE fl273QT> for k = 1, I = 0, with n = 1. 
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FIGURE 4. Rescaled solution Y(y) of the ODE (jlOOj) for k = 1, I = 0, with n = 0.7. 



oscillatory linear eigenfunction of the operator B in (j!.12p 

Mv) = 

For k > 1, the solutions Y(y) are oscillatory (changing sign) as y — > y^, so that 
expansions such as H 2 . 3 1 1) must include oscillatory components as an extra multiplier: 



(2.38) 



Y(y) = (V - yo) a ((f>(s) + o(l)), where s = ln(j/ - y ), 
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FIGURE 5. Rescaled solution Y(y) of the ODE ([235]) for k = 1, I = 1, with n = 3. 
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FIGURE 6. Rescaled solution Y(y) of the ODE ([235]) for k = 1, I = 1, with n = 4. 

and 0(s) is a periodic solution of a (2fc + l)th-order ODE. Examples of such oscillatory 
structures ( 12.381) for k = 2 are presented in [TH pp. 186-192] for k = 2 and /c = 3. 

We will use a similar asymptotic approach in Section [331 in the opposite limit y — > +oo. 

3. On some mathematical aspects of similarity profiles 

3.1. Local existence and uniqueness for k = 1. We need to show that the above 
numerical construction can be justified, by carefully verifying some rigorous aspects of the 
asymptotic analysis. We first apply a fixed point approach for the equivalent nonlinear 
integral equation to prove the expansion (12.331) near the finite left-hand interface. 
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FIGURE 7. Rescaled solution Y(y) for k = 1, I = 2 of the ODE (pT3T|) with n = 3. 



To this end, we look at our integrated second-order equation (I2.30P for k = 1, assuming 
that Y(y) > is strictly monotone increasing (and hence non-oscillatory) sufficiently close 
to the interface at y = y$ < 0: 



J_ y™+iy f or y > o. 



(3.1) y" = — l-lyl-s+ryj, or Y" = 

v ' n+3 1 1 a n+3 

Therefore, we can rewrite our derivatives of Y(y) in terms of the inverse function y(Y): 

Y ' = % = Tryy and Y " = ^(7) = 

So, now we can reduce our differential equation ( 13.1 ft to its equivalent integral form, 



-Jo Yn+iy 
n+3 



2(S/') 2 



(3.2) 



n+3 



y(Y)= y + J , 

\l 2fy(s)s^TTds 



-^fy(s)s«+i ds 



dr = M(y). 



We next prove the following property of the integral operator M in (13. 2ft : 

Proposition 3.1. For small 5 > 0, M(y) is a contraction in C[0,S], with the sup-norm, 
and therefore admits a unique fixed point y{Y) > on (0, 5) giving the unique positive 
solution of the ODE ( \3.1\) on (yo, yo + e) with some sufficiently small e = e(5) > 0. 

Proof. We need to show that, for M(y) to be a contraction, in the C-metric, 
(3.3) ||M(C 2 )-M(Ci)|| </i||C 2 -Ci||, 

for some constant /i = /i(<5) G (0, 1). It is easy to see that M : Z$ — > Z$, for the space Z$ 
of continuous functions, Z$ = {C00 G C[0, 8], ((0) = y }, with the sup-norm: 

= sup ye(M) \£{Y)\. 

14 



Now, take arbitrary d(Y),^ 2 (Y) e Z$. Then from (13.21) we have that 

||M(C 2 ) - M(Ci) || = J^f\\ ( / ^)~" - ( / Ci(s)s^ da)"' || dr, 

v o 

where we use the simplified notation for the integral J, without any limits of integration, 
to mean J Q r . This equality can now be written as 



||M(C 2 )-M(Ci 
Denoting the exponent v = ^j-j-, we see that 
||M(C)-M(Ci)||< 



(/C2(s)s^ds) 2 -(/Cl(s)s^FTds) 
(/C 2 (s)s^ds) ^ + (/Ci(s)sWT ds) 



dr. 



71+3 

2 



n+1 
Y 
I 



f[C2(s)-Cl(s)]s^ ds 



dr. 



Since we deal with sufficiently small values of Y, so that always Ci,2(s) ~ Ho, it is easy to 
estimate in the denominator to get that 



Y 



||M(C a )-M(Ci)||<A*o/ 



IIC2-C1II Jfl"d 8 



Y 



< M0IIC2 — Cill / 





n + 2 



/ s" ds J s" ds[( / a v ds) 7 + (/s"ds) 7 ] 



dr 



Y 



n+2 



n+2 n+r_ n+» dr-< M0IIC2 - Cill Jr- 2 ( „+i) dr < yUo || C2 _ Cl ||y2(n+ 1) 

Q 7- n+1 |. n + I r 2(n + l) q 



Here, /x is a constant dependent on n and y - Since we take | |^/o | < V < 1 2/0 1 ■» then we 
have that \Y\ < 1, and so, fixing Y G [0, 1q], with \i = /x |^o| < 1> we have that (13.31) 
holds true. Hence by Banach's Fixed Point Theorem [TJ p. 39], M(y) has a unique fixed 
point in Zg. □ 

Similarly, we prove local existence and uniqueness of a positive solution for I — 1, 2. 

For k > 2, the solutions are oscillatory (changing sign) close to interfaces, so that a 
contraction approach does not apply, and, as we have mentioned, we need to use other 
techniques of asymptotic analysis in both limits y — >■ y$ and y — >■ +00 (the latter one 
includes k — 1). Such oscillatory structures near interfaces have been thoroughly studied 
for higher-order thin film equations; see [2J [3]. Examples of such oscillatory patterns for 
various NDEs can be found in [T?l Ch. 4], so we do not address these questions here 
anymore, and concentrate on another and more difficult limit. 

3.2. Global existence and uniqueness for k — 1. We restrict our attention to first 
nonlinear eigenfunction Y Q and prove the following: 



Theorem 3.1. For any n > and a fixed interface point y < 0, the problem (I2.30|) . 
(I2.33P admits a unique solutions Y (y), which is infinitely oscillatory as y — » +00. 



Proof. Once the local existence and uniqueness have been established in Proposition 13.11 
its unique existence on the whole interval (yo, +00) follows from elementary checked local 
extension properties for the ODE (I2.30p . which is shown not to admit strong singularities 
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("blow-ups") at any finite point. Oscillatory character of such a solution will be shown 
below. □ 

3.3. Local and global behaviour of nonlinear eigenfunctions as y — > +00. At this 
moment, we do not know the behaviour of nonlinear eigenfunctions fi(y) (or Yi(y)) for 
1/ > 1. Since the first eigenfunction Yo(y) is expected to converge to our linear kernel 
F(y) as n — > + (the Airy function (I2.35P for k — 1), we also expect to have reasonably 
similar behaviour as y — > +00. Recall that the rescaled kernel F(y) = ipoiv) nas the 
following oscillatory slow algebraic decay as y — > +00 [U § 2.2]: for some constant cGK, 

(3.4) F(y) ~ y w cos (d k y^ r + c) as y -)• +00, where d k = 2k (^n) 2k ■ 

Further eigenfunctions ipi{y), given by the derivatives (I2.13p . have the corresponding 
asymptotics via differentiating d3.4j) . so these are much more oscillatory and unbounded 
for y ^> 1. On the other hand, the eigenvalue problem (12. lip admits polynomial solutions 
(cf. f T2TT5|) for B*) 

(3.5) My)~y l + ... for \i = ^ 1 = 0,1,2,.... 

Since ipi ^ L 2 p , these are not proper eigenfunctions. However, this shows that the equation 
(12.111) admits polynomially growing non-oscillatory solutions as y — > +00. 

Due to the complicated nature of the nonlinear equations for n > 0, it is not that 
easy to predict possible behaviour of solutions as y — > +00, where we expect them to 
be oscillatory as in the linear case n = 0. We first study existence and nonexistence of 
non-oscillatory solutions, which mimic the polynomials (I3.5p . Without loss of generality, 
we consider the ODE (I2.2ip for the first nonlinear eigenfunction Y Q (y). 

Proposition 3.2. For any n > 0, the ODE (I2.2ip for even k = 2,4, ... admits a bundle 
of algebraically growing solutions as y — > +00 

(3.6) Y(y) ~ ±Ay m , where m = 2k + gf and A = A(n, k) ^ 0, 

and does not admit such non-oscillatory solutions for odd k — 1, 3, ... . 

Proof. As a formal calculus, we substitute into (12.211) the asymptotic expression from 
(13. 6p to get by balancing the leading terms: 

(3.7) \A\~^i = {-l) k m{m - l)...(m -2k + \)[{2k + 1) + n], 

whence the explicit expression for A and the result. □ 

Thus, for odd k, all the admitted behaviour as y — )■ +00 are not of algebraic growth, 
and, more plausibly, are oscillatory (as well as for all even k). To show this, we separately 
consider the particularly interesting case k — 1. 

Proposition 3.3. All the solutions of the ODE ( \2.21\i for k = 1 are oscillatory as y — > 

+oo, i.e., have infinitely many sign changes in any neighbourhood of +00. 
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Proof. Assume first that, for k — 1, 

Y{y) ~^ +°° as V ~ * +00. 

Then (TCTj) yields 

(3-8) Y" = -^yY^ «-^ 3 for y> i y (y) « 

whence the contradiction. Further cases are studied similarly. The oscillatory character 
of all the solutions follows form the ODE in (I3.8p . since signF" = — signF (cf. Y" = —Y 
for harmonic oscillations). □ 

3.4. A priori bounds for Y(y) for k — 1: nonlinear oscillatory tail. We continue to 
study oscillatory properties of nonlinear eigenfunctions Yi(y), as y — > +00, for the basic 
case k — 1. Let us fix two successive local extremum points 1 <C y\ < y 2 of Y(y), where 
Y'{yi) = Y'{y2) = 0. Let us characterize the size of oscillations of Y(y) at those points. 

We initially look at the equation for the first nonlinear eigenfunction, with the eigen- 
value ao(n), for k — 1. Multiplying ( 12 .30 j) by Y' and integrating over (2/1, 2/2) yields 

m Y "Y' = -4, P 2 \Y\-^iYY'y. 

Simplifying this, we see that 

1 rcy'iT 2 = (n+1) r m (\Y\^Yv 

2 IV J Jj/i (n+2)(n+3) Jj/ivl 1 I J 



After integrating by parts, we obtain 

1 \( Y ') 2 ] V2 = (ra+1) \\Y\^v] y2 + (n+1) P 2 |y|Sr 

2 IV I \y x (n+2)(n+3) U 1 I y J j/i (n+2)(n+3) J J/i 1 1 I 

However, since we are looking at extremum points, where Y'(yx) = Y'^y?) = 0, 
[{Y'f]l = and \Y{y 2 )\m y2 - \Y{y x )\m yx =J%\Y\& dy. 

Since /^lyi™^ 1 dy > 0, we have 

71-1-2 TL I 2 I V/" / \ I TL-\-2. 

(3.9) |y(y 2 )|^y 2 > |F(yi)|~yi, or, on rearranging, (jyjjgf) n+1 > (g). 

This gives a lower estimate on the character of oscillations. 

Let us now look at our second nonlinear eigenfunction with a x (n), where our equation 
is given by (12.241) . k — 1. Multiplying by Y' and integrating yield 

Y'Y"y - {Y'f + ^_|y|-^riyyy = o 

=► |[(y) 2 ]'y - (y) 2 + (^Wy (M^)V = o. 

Integrating between y x and ?/2 again, we have that 

2l\ I J y\yi 2 J 2/1 V 2 J ^ (n+2)(3+2n) LI 1 I # Jyi (n+2)(3+2n) J j/i I 2 I V U 

2(n+l) f2/2|V|5±f, 



_3 m/y/\2 2(n+l) poiy.^ 

2Jt/i^ / (n+2)(3+2n) J 2/1 1 I y 



(I^I^I - |v(yi)l^y?) = 0. 



(n+2)(3+2n) 
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One can see that 

Hi-fi") 2 + jMk O y ^y > °- 

since y > 0. Hence, we have a similar estimate: 

(3.10) \Y(y 2 )\^y 2 2 >\Y( yi )\m y l or (fjgfl)^ > ft) 2 - 

Finally, for the third eigenfunction with c^W; using (I2.26p . k — 1 and multiplying by 
V yields 

Fry - 2(r') 2 i/ + 2ff' + 3^y~^yy'y 3 = o 

=► \[{Y'f\'y 2 - 2(n 2 + <y*y + (w+ ff +3w) (iri^)y = o. 

Integrating over (2/1,2/2) yields 

^[(nviK - /sen 2 * - 2^<n 2 + + j^im im^vX 

It then follows that 

(3-11) (^)^>(S) 3 > and, for any/, (Pffjj) ^ > (g) m 

Thus, the three estimates (13. 9p . (I3.10p . and (13. lip characterize the behaviour of the 
nonlinear oscillatory tail for k = 1, which can be compared with numerical evidence in 
Section 12.51 

3.5. More on oscillatory structure and periodicity. As mentioned before, at least 
for small n > 0, we expect that nonlinear eigenf unctions exhibit, as y — > +00, a behaviour 
which is structurally similar to the linear kernel and its derivatives. We therefore expect 
to have a special oscillatory behaviour for y 3> 1, as we have already seen before. Hence, 
we look to describe this oscillatory structure and, as a first natural attempt, we will try to 
find if these oscillations are given by periodic functions. Let us introduce the oscillatory 
component 0(s) such that, as y — > +00, 

(3.12) Y(y) = 2/ 7 0(s), where s = \iay, 

for some power 7 6 R. Here the term y 1 gives the rate of any growth/decay of the 

2k— 1 

oscillations and may be compared to the controlling factor y ~ found in the linear 
asymptotics (13.41) for n = 0. 

Let us begin with the simpler case k — 1, where substituting into f )3.ip . we obtain 

(3.13) (2/W + ^3 2/ 1+ ^|0|-^0 = 0. 
Expanding this expression and equating powers of y, we find that 

(3-14) 7-2 = 1 + ^ =► 7 = ^- 

This gives us a second-order ODE for the oscillatory component: 

(3.15) P 2 (0)=0" + (2 7 -l)0' + 7(7-l)0 = -^|0r^0 in R. 
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A most typical orbit describing oscillations should be a periodic orbit of (I3.15p . However, 
since 7 > in (I3.14p . the behaviour (13.121) describes unbounded oscillations as y — > +00, 
which are not acceptable for a source-type solution; cf. the linear decaying one (13. 4p for 
n = 0. Hence, we conclude that oscillatory behaviour as y — > +00 is not given by periodic 
oscillatory components <f)(\ny) as in (I3.12p . 

Similarly, we arrive at a contradiction, applying (I3.12p to the general equation ( I2.29[) : 

{-ifDf +\ y ^) = {2k+1) \ {l+l)n y{y^\-^)' + (^TiftW V^M-^- 
Balancing the polynomial terms yields 

(3.16) 7 _ ( 2 * + l) = -J_ => 7 = (2k+ T +1) > 0- 

This leaves us with an ODE of the order 2k + 1, for 0(s): 

(3-i7) (-i) fc iw</>) = (2fc+1) i (i+ i)n + 1 i<r^- 

Here, P2fc+i(0) is a polynomial operator on 0, induced by the term P^ fe+1 (y 7 0). For the 
case k = 1, the polynomial operator is given by 

P 3 (0) = 0"' + 3 ( 2 "+ 3 ) (f)" _|_ 9n 2 +7n+27 0/ + 3(»+l)(2n+3)(n+3) ^ 

In general, P2fc+i(0) is defined by the recursion 

iW(0) = S p 2fc-i(0) + (7 - 2A;)^P 2fc _ 1 (0) + (7 - 2k + l)£P 2fc -i(0) 

+ (7-2fc)( 7 -2A; + l)P 2fc _ 1 (0). 

As usual, periodic solutions of (I3.17P are at most simple oscillatory components. How- 
ever, as before, since 7 > in ( 13. 1 6j) . oscillatory structures of the form ( 13 . 1 2j) . for any 
k > 2, are not applicable, at least for the first nonlinear eigenfunction Y (y), which is 
assumed to be integrable as y — > +00. 

For higher-order nonlinear eigenf unctions Yi{y), with / > 1, proving existence of periodic 
solutions of the ODE (I3.17P is the first step in understanding the oscillatory behaviour. 
This is a difficult mathematical problem, which nevertheless can be solved for orders of 
k that are not too large. We refer to [21 [3], [Tl] for key references and recent results on 
existence- uniqueness of periodic solutions of even-order ODEs such as (13.171) . which occur 
in parabolic thin film theory. We also refer to § 4] for existence results of periodic 
orbits for oscillatory solutions for the PME-4 (I2.ip . m = 2. 

Overall, we rule out the "periodic" structures (13.121) as y — > +00 for the first nonlinear 
eigenfunction Y (y), which is expected to have an oscillatory decay and be integrable 
(not in the absolute sense, i.e., it is not measurable there). Therefore, in this case, 
the oscillatory behaviour may be more complicated and corresponds to not that easy 
"nonlinear focus", which we are going to catch using numerical methods. For I > 1, such 
a behaviour (13.121) with almost periodic oscillatory components <f)(lny) is still plausible 
(but remains rather suspicious). 
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3.6. Branching of nonlinear eigenfunctions at n — 0. As we have promised, we now 
apply another classic idea to trace out the behaviour of all the nonlinear eigenfunctions 
for small n > 0. Namely, we are going to show that there exists branching at n = + 
of solutions with respect to the parameter n. In other words, we show that, as n — > 0, 
there exists certain convergence to solutions (driven by the eigenfunctions of the linear 
operator B in QUO) ) of the the LDE (fT~13]) . 

To this end, let us look at the general ODE given by f 12 . 5 j) . We first expand \f\ n to 
formally get 

(3.18) \f\ n = 1 +nln|/| +0(n 2 ). 

This is pointwise and uniformly true in any bounded positivity subset {|/| > 5o > 0}. 
However, we are not at this moment going to discuss a rigorous functional meaning of 
this expansion for changing sign functions f(y) defined in the whole R. Note that (13. 18)) 
can then be understood in a weak sense, which may be sufficed for passing to the limit in 
the equivalent integral equations; see (21 § 7.6] for asymptotic details. 
Thus, using the formal expansion (I3.18p . ( 12. 5 j) reduces to 

(_ 1 )fc+l jD 2 fe+1 [ (1 + n m IJDJ] + fy + af + 0(n 2) = o. 

Expanding coefficients for small n > yields 

(3.19) (B - \J)f + {-lf^Df+\n In |/|/) + (a - Jg.) / - ^ fy + 0(n 2 ) = 0, 

where B is the linear operator (11.121) and = — 2fc+T ^ s ^ s (^ + -^)^ n eigenvalue. 
For I < 2k + 1, we can find our eigenvalues a/(n) explicitly as in (12.27)) . so that 

a,(n) = m = J±L\l + n (' +1 )1 _1 = h - n (' +1 )1 + 0(n 2 ) 

UL l\"'J (2fc+l)+n(/+l) 2fe+l L 1 ^ 2fc+l J 2fc+l l L 2is+l J T ^ ^ 

Then (13. 19)) reduces to 

(B - \il)f + (-l)^Df+\n In |/|/) - ^ f - fy + 0(n 2 ) = 0. 

Hence using the Lyapunov-Schmidt method [23] by setting 

(3.20) / = ^ + W>, + 0(n 2 ), 

we obtain, within the order 0(n), the following inhomogeneous equation: 

(3.21) (B - A,/M = (-^^^(lnl^iVz) + + j^^iV = h. 

Using Hermitian spectral theory for B and completeness-closure of the eigenfunctions 
subset $ = {i[)i}i>o |H § 4], for the unique solvability of ( 13.211) for <pi, it now remains to 
demand that the right-hand side h is orthogonal to i.e., 

(3.22) (Mz% = 0. 

Here we have to use the corresponding indefinite metric (•, •)*, in which the pair {B, B*} 
comprises the operator B and its adjoint (this metric can be reduced to the standard 
dual L 2 -one, [U § 5]). Here, (13.22)) is known as a scalar bifurcation equation in the classic 
Lyapunov-Schmidt method [23J. 
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We then use the adjoint polynomial eigenf unctions ip* given by (|2.15p . Then, for I < 
2k + 1, we have that the generalized Hermite polynomials are simple [U § 5], 

r l {y) = ^y l (0<Z<2fc+l). 
Hence, for all I <2k + l, (13T22]) is indeed valid: 

<*> = 7vJ [(-l) fe ^ fe+1 dn + (-1)' W + (-1)' #/ m ] dy 

= 7nJ [(-l) fc ^ +1 (ln MlnW + dy = 0. 

Recall that, for I > 2k + 1, we do not know nonlinear eigenvalues az(ra) explicitly. In 
this case, we expand au(n) as follows: 

(3.23) Oiiin) = a + a.\n + 0(n 2 ), 

where «o — cty(0) i n (I2.12p comes from linear Hermitian theory, and «i is a new unknown. 
As before, we use (13. 19[) and now we substitute (13. 23 p . as well as ( I3.20p . to obtain 

n(B - \il)<f)i =(-1)^ + 1 (nln |^|^) + ("o + na x - £L) ^ 

+ "(«o - St) Vfo + °(™ 2 ) = °- 

Equating as usual the terms of the order 0(n), we can find the value of ao, with 



- St = 0, and n(B - Xil)4>, 



(3 24) 

= {-ifDf +1 (nln |^|^) + na^ + n(a - ^L) 0i - ^ 

Hence, substituting into (I3.24p and passing to the limit n — > + yield 

(B - = {-l) k Df + 1 (ln + - (ifelF My = h - 

Then, the orthogonality condition (I3.22p becomes an algebraic equation for a\ in ( I3.23p . 
Namely, taking the inner product with ^* and noting that {h, ip^)* = 0, (ipi, ■0*)* = 1 yield 

(3.25) a x = -{{-l) k Df+\\n\^{) - 

Thus, this uniquely defines the second coefficient oti in the expansion (I3.23P and then, 
as usual, (I3.25P gives a unique function ipi in the eigenfunction expansion (I3.20p . etc. 

In the analytic or even finite regularity cases, solvability conditions and existence of 
expansions such as (I3.23H usually rigorously justify the actual presence of branching. 
Our case is more delicate in view of the "weakness" of the expansion (I3.18p . However, 
for the variable Y = \f\ n f, the expansion (I3.18P is easier to justify, especially now the 
equation becomes semilinear and can be reduced to an integral equation with compact 
Hammerstein-Uryson-type operators. Therefore, in the present case, a rather full justifi- 
cation of the n-branching method, though being rather technical, is doable and does not 
represent a principally non-solvable problem of nonlinear integral operator theory. 
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4. The nonlinear limit n ->• oo: = 1 and Z = 



While we dealt before with a "homotopy path" construction of nonlinear eigenfunctions 
in the limit of small n — > + , we now consider the opposite "highly nonlinear" limit 
n — > +oo, which also helps to understand properties of the nonlinear eigenvalue problem. 



4.1. Reducing to an algebraic problem. As we have seen in Section [3.6[ as n — > 
0, there appears a direct connection of all the nonlinear eigenfunctions with the linear 
ones associated with the LDE (11.1 3ft . Now, we are going to perform the opposite highly 
nonlinear limit n — > +oo for nonlinear eigenfunctions of the NDE. Rather surprisingly, 
it turns out that this "limit nonlinear case" admits a more profound analysis, and for 
Z = 0, 1, 2 we are able to tackle the nonlinear eigenvalue problems by a simpler geometric- 
algebraic, allowing us to obtain a number of analytical and explicit expressions. 

Consider the ODE (|2.21|) . with k — 1 and where Z = 0, i.e., (I2.30p . Since we are dealing 
with n ^> 1, it is necessary to scale out any coefficients containing large n's. In order to 
do this, we let 

(4.1) Y(y) = CY(y), where C = C{n) > is a constant. 
Substituting g3J into ([2T30]) yields: 

(4.2) CY" = \C\-^C\Y\-^Yy =► C = (n + 3)-^, 
so that, on scaling out such a C(n), we obtain the rescaled ODE 

(4.3) Y" = -\Y\-**Yy. 

Here we can pass to the limit n — > +oo, since the only n-dependent exponent satisfies 
— — > — 1. At n = +oo, the ODE becomes simpler and contains a bounded discontin- 
uous nonlinearity: 



(4.4) B DO (y) = Y" + sign Yy = in 



Of course, passing to the limit to arrive at (14.41) might be a delicate mathematical problem. 
A potentially dangerous situation occurs in those subsets, where Y vanishes. However, if 
both Y(y) and Y(y) have a. a. zeros transversal in the natural sense, then passage to the 
limit is straightforward. A sufficient "transversality" of zeros of the limit function Y(y) 
can be checked a posteriori, after completing our algebraic construction. 

Solving "almost linear" ODE (14.41) . we find expressions dependent on the sign of Y: 



(4.5) 




Y +(y) = -ly 3 + ciy + c 2 , 

Y-{y) = \y z + d iy + d 2 . 



Here c\, c 2 , d\, d 2 are all constants, not necessarily positive ones. Knowing the conditions 
of continuity for the function Y(y) and Y'(y), we must have that all one-sided limits 
coincide, i.e., 

(4.6) Y + = Y- and Y' Jr = Y'_, at any zero, where Y = 0. 
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Let the points {y = yt}i>o be successive zeros, i.e., Y(yi) = 0, for i = 0,1,2,.... 
Hence, for Y + (yi) given in fj4.5j) . for a fixed isolated zero with an i > 1 (as usual, y < 
corresponds to the left-hand interface), 

Y+iyi) = -\ yf + cum + = o. 

We now rearrange this to find one of the unknown parameters, in terms of and the 
parameter c 2 j, such that cu = | yf — We also have that Y+(yi) = Yl(j/ 4 ), hence 

-|y 4 2 + c u = |y 2 + dii. 

From this, we find a second parameter in terms of yi and c 2 i, where du = — ^yf — ^l. 
Now, we see that from Y + — Y_, 

ci»2/» + diij/i + c 2i + d 2 j = 0. 

So, substituting in known values, we have our third parameter d 2i given by d 2 i = ^yf + c 2i . 
We now see that, after substituting values for cu, du, and d 2i , (|4.5p can be written as 



Y > : y + (y) = -i y 3 + (± % 2 - 2* )y + c 2i , 
y < : Y_(y) = W- (lyf + ^)y + |yf + c 2j . 



(4.7) 

11 .-"II- V I ;; ) — - - I - ... ■ 

2/i /tf 3 

From the above, it is noted that yi ^ 0, for any z, unless c 2 j = 0. 

4.2. Existence, uniqueness, and zero properties of Y (y). We now resolve the al- 
gebraic system to get a rather complete description of some important properties of the 
solution Y (y) obtained by such a "geometric approach. Recall first the scaling invariance 
of the ODE (03}: 

(4.8) Y(y) is a solution ==>- ±a 3 Y(^) is a solution for any a > 0. 
Therefore, choosing the interface at yo = — 1, we put two conditions there 

(4.9) Y(-l) = F'(-l) = 0, 
and prove the following: 

Theorem 4.1. The problem ( I4.4p . (14. 9p admits a unique nontrivial solution Y (y), which 
has transversal zeros {yi}i>i such that 

y = -1, yi = 2\y \=2, y 2 = 3^2 - 1 = 3.2426... , and 
J/i+i = 2 for any i > 2. 

Proof. As we have promised, we construct such a solution using pure algebraic manipula- 
tions. Let us begin with the first interval of positivity (yo = —l,yi), where according to 
f)4.9p . the solution reads 

(4-11) Y + {y) = -\{y + lf(y- yi ). 
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Since (14. 5p implies no quadratic term ~ y 2 in the cubic polynomial, this uniquely gives 
Hi = —2y = 2, and hence 

(4.12) Y + (y) = -±(y+iy(y-2) >0 on (-1,2). 

On the next interval (y\ = 2,y 2 ), the negative solution takes the form: 

(4.13) Y_(y) = l(y-2)(y-y 2 )(y + Cl ) < 0. 
Similarly to the above, we then conclude that 

(4.14) Cl =2 + y 2 . 

Then matching of the first derivatives (14.61) at y — y\ — 2 implies 

(4.15) -|3 2 = |(2- 2/2 )(4 + 2/2 ) =^ yf + 2^-17 = ==► ^ = 3^2-1, 

and this procedure can be continued. 

Consider an arbitrary interval (yi-i,yi) of positivity (or negativity), i > 2, with 

(4.16) Y+(y) = (y-yi-i)(y-yi)(-ly + c i - 1 ), where c i -i = ~(y i +y i - 1 ). 
Similarly, on the next negativity (or resp. positivity) interval (yi,yi + ±), 

(4.17) Y-(y) = (y-Vi)(y-yi+i)(^V + Ci), where c % — \ (Vi + Vi+i)- 
Therefore, the matching at y — y\ yields 

(Vi- Vi-i) ( - | Vi - I (.Vi + Vi-i)) = (Vi - Vi+i) (lVi + 1 (Vi + Vi+i)) 

(4.1oJ 

y*+i + j/ti/t+i - %i + ViVi-i + Vi-i = o, 

whence the final result in (14.101) . Thus, the unique solution can be extended indefinitely 
for arbitrary y ^> 1 and is infinitely oscillatory as y — > +oo. □ 

Remark. The quadratic equation in ( 14. 18ft reduces to a 2D linear discrete equation: 

(4.19) Oij = yiyj =S> Oi+M+i + ~ + a M-i + "i-i.t-i = for z > 2. 
It is not difficult to find some particular solutions: 

(4.20) a iJ = ^ ^ + V-V + /^ + l y±M = 2M± X"+iP - 

Therefore, denoting by .M a proper subset of parameters /i such that {fi l }^eM is com- 
plete/closed, the general solution of (I4.20p is represented by a converging infinite series 

(4-21) aiJ = £ M6 M<V^ 

where {C M } are constants and i/ = z/(/i) take values according to (14.20)) . Here, (14.21)) 
is a discrete analogy of eigenfunction expansions of solutions of a linear PDE with two 
independent variables (x,t). A proper posing "boundary conditions" for ( 14. 19)) to specify 
the corresponding Sturm-Liouville problem and next initial conditions to get the corre- 
sponding eigenfunction expansion (I4.21j) is a difficult and uncertain proble which will 
unlikely provide us with any useful finite explicit formulae. 



1 It seems, a suitable behaviour at infinity, as i,j — > +oo, of a^j might include something like the 
''minimality" condition in (j2.6p . which is hard to take into account. 
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However, the analytic relationships such as (I4.10p can promise to get extra asymptotic 
properties of the first rescaled eigenfunction Yo(y), especially the decay rate of the minimal 
behaviour indicated in (12. 6p . However, in the present case I = 0, unlike the simpler one 
I = 2 in Section [6j some computations are not expected to be easy all the way, since such 
algebraic relations are quadratic and hence not always explicitly solvable. 

4.3. Numerics for n ^> 1: k = 1 and 1 = 0. As usual, very sharp proper numerics can 
help to detect further properties of Y (y), and hence avoid trying to get too complicated 
and exhaustive results concerning the corresponding algebraic system. Moreover, which 
is even more important, we can also check the character of convergence of solutions as 
n — > +oo, which, after scaling (14. ip . turns out to be rather fast. To deal with the ODE 
f)4.4p . it is possible to just use a simple shooting method using the ODE solver ode45, to 
find suitable profiles and nonlinear eigenfunctions. 

We use a similar shooting method as that applied for the general case of n > 0, set out 
in Section |2~51 We recall that, close to the interface at some y — yo < 0, we look for small 
solutions of Y{y) such that, as y — )■ y$, 

Y(y) = C (y-y )l(l + o(l)), where C = \ \y Q \ > 0, 

and Y'(y) = \y \{y-y ){l + o(l)) (k = 1). 

The proof of this expansion is similar and even simpler than that of Proposition 13.11 

In Figure El we show the general view of the first eigenfunction Y (y) on the large 
interval [yo = —1, 100]. The envelope of the decaying oscillations are governed by the 
algebraic curve 

(4.23) L (y) « ±0.7y~3 as y ->■ +oo. 

It seems that, definitely, this decay can be seen from the above algebraic system. We 
have checked that, for n = 100, the corresponding nonlinear eigenfucntion (after scaling) 
is practically indistinguishable. 

In Figure El we show Yo(y) on a smaller interval y G [—1,7], indicating all the first 
zeros, which will coincide with those given by explicit algebraic expressions in f)4.10p . 

All the computations have been performed with the enhanced accuracy and tolerances 
~ 10~ 10 , so these are quite reliable. Let us present first 15 zeros of Y (y): 

y = -1, yi = 2\y \ = 2, y 3 = 3.2426... , y A = 5.0777... , y 5 = 5.8426... , 
y 6 = 6.5459... ,y 7 = 7.2021... , y 8 = 7.8207... , y 9 = 8.4083... , y w = 8.9697... , 
yn = 9.5086... , y 12 = 10.0279... , y 13 = 10.5299... , y u = 10.0164... , ... . 

4.4. Branching at n = +oo. Introducing the small parameter in (14. 3p for n ^> 1, 

e = 1 — -> + as n -> +oo, 

and performing a standard (formal, as usual, at least in the pointwise sense) linearization 
yield the following problem: 

(4.24) \Y\-^h = |y|*-i = i+einm+o(s 2 ) =^ Boo{Y) — —sigaYy hi \Y\ + 0(e), 
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The first nonlinear eigenfunction Y Q (y), n=°°: oscillations, decay, envelopes 




100 



Figure 8. The first nonlinear eigenfunction Y Q (y) as a solution of the prob- 
lem (@3D, (EL22D; k = l, * = o. 



where is the unperturbed operator in (j4.4j) . Next, studying the branching at e = 
from the first nonlinear eigenfunction V^, we obtain the corresponding linear problem: 

(4.25) Y = Y o + e0 + O(e 2 ) => B'JW = h = -(signY )y In \Y \, 
where the symmetric linearized operator B'^Yq) is given by 

(4.26) B'^Yo) = D 2 y + (signF )V, (signfo)' = S(y - y ) + 2 £ (-l)^(y - y<). 

(i>l) 



Here {yi} are zeros of Yo(y) as explained in Theorem I4.ll 

At this stage, one then needs proper spectral theory for a self-adjoint extension of the 
operator B'^Yq). Then the branching condition reads as the orthogonality 

(4.27) h±herB' oo (y ). 

However, developing such a proper spectral theory faces some hard algebraic difficulties. 
Indeed, looking for eigenf unctions (pk, 

(4.28) ipl - \ k ip k = C k = const. = - (y ip k (y ) + 2 ^(-l)*^^^)) , <fk(yo) = 0, 
for Xk < 0, we obtain the solution and an algebraic equation for such eigenvalues: 



(4.29) <p k (y) =sm{ y /\X k ~\{y-y )) 



^k- E (-l) i 2/iSin(vlAfe|(yi -yo)) = o. 

(i>l) 
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The first nonlinear eigenfunction Y Q (y), n=oo: first zeros 




y 

Figure 9. The first nonlinear eigenfunction Y Q (y) as a solution of the prob- 
lem (S3), ffl~22l) : k = 1,1 = 0. 

According to (I4.29p . we require that each eigenfunction <pk(y) should be purely oscillatory 
as y — > +00 about zero, meaning a "minimal" ("zero-average") behaviour or zero value 
in a weak sense. For A = 0, we take <fo(y) — y — yo (again up to a normalization 
multiplier), leading to another algebraic problem, which is expected to be non-proper since 
the eigenfunctions is not oscillatory and not satisfying the "zero condition" at infinity. 

Proving that the algebraic equation in (14.291) has a discrete family of solutions {A& < 
0} is a difficult open problem. Nevertheless, at least, this analysis shows a principal 
possibility to detect branching of nonlinear eigenfunctions at n = +00, where a simpler 
algebraic treatment is available, and more practical asymptotic and other properties of 
the nonlinear eigenfunctions than in Theorem 13.11 for any n > 0, are known. 

4.5. The higher-order case k > 2 for I = 0. For the general case of the higher-order 
ODEs (12.211) . again for / = 0, the equations as n — > 00 are much the same. Here the 
scaling constant C(n), for Y = C(n)Y is given by 

C{n) = [(2Jfe+l)+n] n ■ 
After scaling, this yields the ODE 

(4.30) B^ k \Y) = {-\) k+1 DfY + signFy = 0, 

which deserves further study by deriving the corresponding algebraic structures. Note 
that, close to the interface at y = y^, for k > 2 (the case k = 1 is not oscillatory, as we 
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have seen in Section [3]), a stabilization to a periodic oscillatory component <j)(s) , with the 
typical structure near the interface at y = 



(4.31) Y(y) = {y -y ) 2k {<f){s) + o(l)), where s = - ln(y - y ) ->■ -oo, 

of solutions of ( 14.30P is most plausible. See [2j 13] , where such oscillatory sign changing 
solutions such as (14.311) . with a periodic component <p(s) satisfying a nonlinear higher- 
order ODE, have been found for the fourth- and sixth-order thin film equations. 

However, explicitly solving the ODE (I4.30p for any k > 2 in the positivity and negativity 
domains leads to much more complicated algebraic systems, which do not admit such a 
clear explicit resolving as for k = 1. Moreover, even shooting numerical methods lead to 
difficult and unclear results, so we do not present such an analysis here. 



5. The second eigenfunction Y\{y) for n = +oo: algebraic approach 

Using the same method, the equations governing the behaviour as n — > oo, relating to 
the second nonlinear eigenfunction Yx(y), can easily be found. Namely and analogously, 
for I = 1 in the ODE (12.241) . with k = 1, the limit equation at n = +oo is given by 

(5.1) Y"y -Y' + sign? y 2 = 0. 

This ODE (15. ip can be written in a singular Sturm-Liouville form 

(5.2) (!)' = _ sign f, 



where the weight p(y) = i ^ L\ oc . Therefore, y = is a singular inner point, where an 
extra condition must be posed. This is done by checking the functional weighted L 2 -space 
corresponding to the operator in (I5.2p and its available asymptotics as y — > 0: 

(5.3) J p(y)(Y'(y)) 2 dy<oc =^ Y'(0)=0. 

Next, integrating (15. 2 j) twice yields, in the positivity and negativity domains, the following: 

(5.4) Y± (y) = T v i + ay 2 + b, a,b,eR. 

Note that, unlike (14. 5 p for 1 = 0, here the linear term ~ y is absent in the cubic polynomial. 
The analysis of the algebraic matching system corresponding to (15.41) is similar in many 
places but a couple of ones, which we will concentrate upon now. 

Theorem 5.1. The problem ( 15. 2ft . ( 14.91) admits a unique nontrivial solution Y\(y) with 
transversal zeros at {yi}i>i, where 

(5.5) y = -1, Vl = \, y 2 = ^ = 2.927050..., 

and triples of zeros {yi-i,yi,yi+i}i>2 satisfy a cubic homogeneous algebraic equation, 

(5.6) yi-iy 2 i+1 + yiy 2 i+1 + yl_ x yi + yf-iUi+i = Vi-iViVi-i + vhi+i + Vi-iVi + vf f or i > 2 - 
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Remark. (15.61) reduces to a 3D linear discrete equation for a^j^ = yiVjUk, but, similar to 
(I4.18p . this does not essentially help to get any finite explicit expressions for zeros {?/«}. 

Proof. The first step is elementary: on (— 1, 2/1), in view of (15. 4 p (cf. also (15. 3p ). 

(5.7) Y + (y) = (y + l) 2 (y- yi )(-l), where y x = \. 
Next, on (yi,y 2 ) by (El, 

(5.8) Y-{y) = (y-y 1 )(y-y 2 )(ly + c 1 ), where Cx = \^ 
so that the matching (14. 6 j) at y = y\ yields 

(5.9) (2/i + l) 2 (-|) = (yi-2/ 2 )(|yi + c 1 ) 4y 2 2 - 10y 2 - 5 = 0, 

whence the desired value of the root y 2 in (15. 5p . 
Similarly, in the general case, on (y^,^), 

(5.10) y_(y) = (y - y t -i)(y - Vi)(\y + C;_i), where q_i = | ff^gr , 
and on (yi,y i+1 ), 

(5.11) r + (y) = (y-y,,)(y-y m )(-|y + Q), where Q = -§ -»±-. 

Then the standard matching via (14.61) of such Kt(y) yields the homogeneous cubic alge- 
braic equation (I5.6P for the zero triple {yi-i, yi, yi+i}- □ 

Figure [10] shows the general structure of the second eigenfunction Y\{y) on the large 
interval [yo = —1,200]. The envelope of the decaying oscillations are governed by the 
algebraic curve 

(5.12) Li(y) w ±0.89 as y -)> +oo, 

which can be associated with the algebraic manipulations involved. 

The next, Figure [TTl shows first zeros of Yi(y) on the interval [—1,5]. 

In general, for all values of k > 1, we have the limit equation at n = +oo for / = 1 in 
the form: 

{-l) k+l Df Yy + {-l) k Df~ 1 Y + signy = 0. 

This can be integrated once, however, one cannot expect any reasonably easy algebraic 
manipulations leading to some explicit representation of asymptotic properties of Y\(y) 
for any k > 2. 

6. The third eigenfunction Y 2 (y) for n = +oo: algebraic approach 

Similarly, for the third eigenfunction, where I = 2, we have, in the lower-order case 
k = 1, the limit ODE 

(6.1) Y"y 2 -2Y'y + 2Y + signY y 3 = 0. 

The corresponding Sturm-Liouville form 

(6-2) (g)' = -Af-I sign f, 
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The second nonlinear eigenfunction Y^y), n=°°: oscillations, slow growth, envelopes 
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Figure 10. The second nonlinear eigenfunction Y(y) as a solution of the 
problem CT) . (14321) : fc = 1, / = 1. 

reveals even more singularity at y = as it used to be in ( 15.2ft for The necessary- 

extra condition at y = is presented below at ( 16.101) . 
Writing (16. ip down for Y±(y) in the form 

(6-3) Y" = T y + 2^y 

and integrating twice yields 

(6.4) Y± (y) = T^ + ay 2 + by, a,beR. 

Surprisingly, the matching at zeros y — y% by using the cubic polynomials (16 ,4p . having 
y = as a fixed zero always leads to a simpler mathematics. 

Theorem 6.1. The problem ( I6.ip . ( 14. 9 p admits a unique nontrivial solution Yi(y) with 
transversal zeros at {yi\i>\, where 

(6.5) 2/o = -1, 2/i = 0, 2/2 = 1 + ^ = 2.4142..., 2/3 = 1 + 3^= 5.2426..., 
and further zeros are given by the second-order linear discrete equation: for any i>3, 

(6.6) 2/i+i - 2 Vi + yi-i = =*► j/i = Ci + C 2 «, w/iere Cj = 1 - 3^2, C 2 = 2^2. 
In particular, the distribution of zeros is uniform: 

(6.7) 2/i+i -V* = 2V2 = 2.8284... for all i > 2. 
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The second nonlinear eigenfunction Y (y), n=oo: first zeros 
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Figure 11. The second nonlinear eigenfunction Yi(y) as a solution of the 
problem CT) . fTC22D : fc = 1, / = 1. 

Proof. On (— 1, yi), in view of (16. 4p . 

(6.8) Y + {y) = y{y + lY{-\) Vl = 0. 
Next, on ( Vl = 0,y 2 ) by (Q, 

(6.9) y_(y) = y(y-y 2 ) (§ y + d), where d = ^. 

Unlike the previous cases, it is not possible to find y 2 and c\ by using the standard 
matching conditions (14. 6p . The point is that the differential operator in (I6.2p is strongly 
singular at y — 0, where the weight p(y) = -j ^ L p (— 1, 1) for any p > 1. 

It then follows from f)6.4p due to the singular setting (I6.2p that the two usual conditions 
at y = such as the values of Y(0) = and of a given "flux" Y'(0) are not sufficient 
to determine a unique local solution for y > and ?/ < 0. To get a unique solution, 
the value of Y"(0) should be prescribed. A further analysis shows that a proper stronger 
continuity condition of matching at y — is necessary and this includes the equality of 
the second-order derivatives: 

(6.10) Y^(0~) = F"(0 + ). 
Overall, this yields the following quadratic equation for y 2 : 

(6.11) F|'(0) = -2 = ri'(0) = 2 Cl - y 2 , c x = ± =^ y\ - 2y 2 - 1 = 0, 
and this uniquely defines y 2 shown in (16. 5p . 
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(6.17) 




Next, we use (I6.9P with the obtained values of y 2 and c\ to match (now, in a standard 
way) with the solution representation on (y 2 ,ys), 

(6.12) Y + (y)=y(y-y 2 )(y-y 3 )(-l), 
to get at y = y 2 

(6.13) 3/2(52/2 + 01) =2/ 2 (2/ 2 -2/ 3 )(-i) =^ 2/3 = 2y 2 + 2 Cl = 1 + 3^2. 
Finally, in the general case, on (yi-i,yi) for i > 3, 

(6.14) y + (y) = - yi _!)(y - - i ), 
and on (yn,y i+1 ), 

(6-15) F_(y)=y(y- 2 / i )( 2 /-y m )i. 

By the standard matching at y = yi via ( 14. 6 p of such Y±(y) yields 

(6.16) y t (yi - y^ x ) ( - 5) = vAVi - Vi+\% 

whence the linear difference relation ( 16. 6ft with some constants C\ and C 2 . These are 
uniquely obtained from the linear algebraic system, 

'y 2 = 1 + y/2= C ± + 2C 2 , 
y 3 = l+3>/2 = C 1 + 3C 2 

completing the proof. □ 

Figure [12] shows the general structure of the third eigenfunction Y 2 (y) on the large 
interval [y = —1,200]. Directly connected with (I6.6p . the envelope of the decaying 
oscillations is linear 

(6.18) L 2 (y) ps ±y as y ->■ +00. 

First zeros of Y 2 (y) on the interval [—1, 7] are shown in Figure [TBI 
For arbitrary k > 2, the limit ODE for Y 2 (y) is 

{-l) k+1 DfYy 2 + 2{-l) k Df^Yy + 2(-l) fc+1 £)J*- 2 y + X y 3 = 0; 

which does not admit a clear geometric-algebraic method of solution. 

Finally, we recall that exotic profiles such as Y (y) in Figure [SJ Yi(y) in Figure ITOl 
and 1^2(2/) in Figure [12] are not just some functions defining some self-similar solutions in 
the nonlinear limit (n = +00) equation (I1.16P for k — 1, but these are the most stable 
asymptotic pattern of such a nonlinear PDE. Indeed, Y (y) is the most stable and is 
expected to attract, as t — > +00, a. a. solutions, excluding only those with data satisfying 
some extra orthogonality conditions. For instance, those having zero mass (or also the 
zero moment for Y 2 (y) to play a role), so that an extra time-scaling is necessary to get 
convergence to the second nonlinear eigenfunction Yi(y), etc. Of course, these stability 
questions are far beyond the scope of this paper, are very difficult, and remain open for 
most higher-order NDEs. 
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7. Nonlinear dispersion equation with absorption 

7.1. A full quasilinear NDE: homotopy deformation to linear PDEs. The final 
natural progression, from the nonlinear model (11. lip , is to go to the odd-order NDE with 
absorption (11.11) . This links up both the nonlinear model and the semilinear one [I], 

(7.1) m = {-l) k+1 Df +l u - \u\ p -\ in KxK + (n = 0). 

However, the NDE (ll.ip is indeed a more difficult quasilinear equation than (17.11) . with 
not that well understood phenomena of "nonlinear" bifurcation and branching. Therefore, 
we should pay here less effort towards establishing rigorously some of the key analytical 
aspects concerning basic similarity solutions of ( II. ip . 

The NDE ( 11.11) . for any n > 0, p> n + 1, after similarity scaling, 

(7.2) =t-5^/(y), where = j^ggLy, 
reduces to the ODE 

(7.3) (~l) fc+1 £>f + <£§£li) f'v+^f- l/r 7 = m K. 

As usual, a proper setting for (17. 3p assumes a finite left-hand interface at some y = y < 
and an admissible oscillatory behaviour as y — » +oo, which was under scrutiny above. 
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The third nonlinear eigenfunction Y 2 (y), n=°°: first zeros 
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Figure 13. The third nonlinear eigenfunction Y 2 (y) as a solution of the 
problem (Oil . (147221) : fc = 1, Z = 2. 

Recall that, using the reflections (12.81) . simultaneously, we construct blow-up solutions 



(7.4) 



u u {x,t) = (T-ty—f{y), y 



where f3 



p-(n+l) 



(T— tp ' K (p-l)(2fc+l) 

(since /3 > for p > n+ 1, formally, this is a single point blow-up) of the NDE with source, 
(7.5) 



(\u\ n u) + (w^w 



so that proper profiles f(y) describe both global and blow-up asymptotics of such NDEs. 

In the case n = in (17. 3ft . we have a simpler semilinear ODE corresponding to the 
model f l7.ip . It is important also to note that, in addition to the homotopy path n — > + 
in ( ll.ip . it is useful to apply an extra limit p — > 1 + . For the case n = and p = 1 in (11. ip . 
we reduce to the linear equation 

Ut = (-l) k+1 Dl k+1 u - u. 



Using the substitution 



u(x, t) = e v(x, t), 



reduces it to the standard linear dispersion equation 

(7.6) v t = {-l) k+l D 2 x k+l v. 

The extra scaling yields a semigroup (a group) with the infinitesimal generator B in f l 1 . 1 2 p : 

(7.7) v(x, t) = t~™+iw(y, t), ?/ = x/t2fe+r, T = h\t =^ w T = Bw. 
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Therefore, linear Hermitian theory developed in [U § 4] (see also [61 §9]) can be used, and 
this leads to an efficient approach based on comparison of linear eigenfunction structures 
of (17.61) and the nonlinear ones of (11.11) . Actually, this implies existence of a certain 
"homotopy" of those PDEs as n — > + (already noted in Section 13. 6ft and as p — > 1 + . 
Overall, this establishes a countable nature of the nonlinear eigenfunction family for (II. ip . 
which is difficult to prove rigorously. We refer to [15] , where such a homotopy is discussed 
for a related fourth-order PDE. 

Thus, the continuous double limit, n — > + and p — > 1 + , by reducing to the LDE (17. 6p 
with a countable set of rescaled eigenfunctions (I2.13p . in general, justifies that the ODE 
(17. 3p admits a countable set of so-called p-bifurcation branches of solutions, which blow- 
up as p — )■ (n + 1) + . Further study of such nonlinear phenomena is of importance. For 
n = 0, this branching phenomenon has been studied in [4]. 

It is worth mentioning that, in general, such a full branching approach assumes existence 
of a countable number of 2D bifurcation surfaces defined in the 2D parameter space 
of (n,p), and the critical point (0,1) is expected to be a complicated singular "cusp" 
induced by such a collection of bifurcation surfaces. This is a very difficult bifurcation- 
branching phenomenon, which in 2D was not fully understood and remains an open 
problem. Therefore, below, we study another formal, but ID, "nonlinear bifurcation" 
phenomenon in order to explain existence of a countable number of p-bifurcation branches 
in the problem (I7.3p . 

7.2. Local "nonlinear" bifurcations at critical exponents pi(n). We consider the 
VSS solutions (17. 2p of (II. ip . and develop a formal nonlinear version of such a ID p- 
bifurcation (branching) analysis for any fixed n > 0. As usual, according to classic 
branching theory [T8~| [23] . a justification (if any) is performed for the equivalent semilinear 
integral equation with compact operators in suitable metrics. For simplicity, we present 
computations in the differential setting, which does not change anything essentially. Note 
that, for nonlinear odd-order operators, some issues of compactness can be rather tricky. 

Let us compare the time-factor structure of the VSS (17.21) and that for the pure NDE in 
(12. 2p . It follows that the critical bifurcation exponents {pi = pi{n)} are then determined 
from the equality of the exponents: 

(7-8) ~iFI = -«' =► ^) = 1 + 4rj 1 = 0,1,2,..., 

where cti(n) are the critical exponents as in (I2.27P obtained explicitly and further fully 
nonlinear ones that cannot be determined dimensionally (via conservation laws). 

In particular, for the semilinear case n = 0, the eigenvalues «;(0) are given by (I2.12p . 
and this leads to the critical bifurcation exponents 

(7.9) W ( ) = 1 + 2M!, Z = ,l,2,..., 

for the semilinear equation (17. ip studied in [JJ § 7.4]. In the present nonlinear case with 
a fixed n > 0, such a standard linearised approach is not suitable. 
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However, first steps of this "nonlinear" bifurcation theory are straightforward. We next 
use an expansion relative to the small parameter e — pi — p, so that, as e — > 0, 

a = = ^TS = ai + ocfe + 0(e 2 ) , 

a _ p-(n+l) ._ \-a x n _ naf n(r 2 \ 
1 ~ (p-l)(2A+l) 2fc+l 2/c+l ^ ^' 

l/r 1 / = l/l* -1 -/ = l/l Prl /(l - e In |/| + 0(5 2 )). 

Note that, unlike the case (I3.18p . the last expansion has a clearer functional validity, 
since at / = 0, there occurs standard issues of convergence, which makes sense suitable 
for passing to the limit in the integral operators. 

Substituting these expansions into (I7.3p and collecting 0(1) and 0(£)-terms yield 

(7.10) A(/, a,) - l/p*- 1 / + eCf + el/l*" 1 /!* \f\ + 0(e 2 ) = 0, C = -gt yD y + afl, 

where A(/, cti) is the nonlinear operator in f 12 . 5 j) . Recall that, at each nonlinear eigenvalue 
a = cti, there exists the corresponding nonlinear eigenfunction f) such that (12. 71) holds. 
At least, we are going to use this conclusion, which was not completely proved. The fact 
is that the operator A(f,a), with a = on in (17. lOf) of the rescaled pure NDE, correctly 
describes the essence of a "nonlinear bifurcation phenomenon" to be revealed. 

To this end, we use the additional invariant scaling of the operator A(/, a) by intro- 
ducing the new unknown function F(-) as follows: 

(7.11) f(y) = bF(y/b»^), 

where b = b{e) > is also an unknown parameter to be determined from a scalar branching 
equation and satisfying 

(7.12) b(e) ^0 as e 0. 

Substituting (17. lip into (17.101) and, under natural (but not proved) regularity assump- 
tions on such expansions, omitting all higher-order terms (including the one with the 
logarithmic multiplier In yield 

(7.13) A{F, ai) - If-^F^F + eCF = 0. 

Finally, we perform linearization about the nonlinear eigenfunction fi(y) by setting 

F = fi + Y. 

This yields the following linear non-homogeneous problem: 

(7.14) A'(f hai )Y = V'-^f^fi-eCfi. 
Here, the derivative is given by 

(7.15) A'(/,,a,)Y = {n + l^-lf^Df+W^Y) + ^Y'y + a l Y. 

As usual in bifurcation theory, the rest of the analysis crucially depends on assumed 
good spectral properties of the linearised operator A'(fi,ai), which are not easy at all, 
and, in fact, are much more complicated than that for the pair {B,B*} in [4J. Many 
aspects of such a theory remain quite obscure. However, we proceed to explain the key 
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final results on possible bifurcations, and follow the same lines. A proper functional 
setting of this operator is more understandable in the present ID case, where, using the 
behaviour of fi(y) — > as y — > y£ and y — > +00, it is, at least formally, possible to check 
whether the resolvent is compact in a suitable weighted L 2 space. In general, this is a 
difficult open problem, especially since we are not aware of precise asymptotic properties 
of all the eigenf unctions fa. 

We assume that such a proper functional setting is available for A'(fa,ai). Therefore, 
we deal with operators having solutions with "minimal" singularities at the boundary 
point of the support at y = y < 0, where the operator is degenerate and singular. The 
same is assumed at y — +00, where the necessary admitted bundle of solutions should be 
identified to pose singular boundary conditions; see Naimark's monographs on ordinary 
differential operators [121 ED] as a guide. 

Namely (cf. [4, § 4]), we assume that the linear odd-order operator A' (fa, an) has 
a discrete spectrum, and a complete and closed set of eigenfunctions denoted again by 
{ip-y}- We also assume that the kernel is finite-dimensional and we are able to determine 
the spectrum, eigenfunctions {ipt}, and the kernel of the "adjoint" operator (A'(fa,ai))* 
defined in a natural way using the topology of the dual space L 2 (or, equivalently and 
possibly, of a space with an indefinite metric) and having the same point spectrum. 
The latter is true for compact operators in suitable spaces in more standard setting, [T7] 
Ch. 4]. We also require that the bi-orthonormal eigenfunction subset {^7} of the operator 
A'(fi,aci) is complete and closed in a suitable weighted L 2 -space (for n = 0, such results 
are available [4]). Note that, often, this "spectral collection" is too exhaustive in nonlinear 
operator theory; see Deimling [TJ p. 412] for most general bifurcation results. 

Thus, by a typical Fredholm-like alternative, the unique solvability of (I7.14p requires 
the orthogonality of the inhomogeneous term therein to the ker A'(fi, a{). For simplicity, 
let it be ID with the eigenfunction <pi, so the right-hand side in (17.141) satisfies 

(7.16) W l - l \fi\ Pl - l fa-eCfa 1 kerA'(./i,a,) = Span{0,}. 

Then, multiplying (17.161) by <\)\ in L 2 (or within the equivalent indefinite metric) yields the 
orthogonality condition (Lyapunov-Schmidt's algebraic branching equation [231 § 27]): 

(7.17) b*- 1 (\far- 1 f h tf)=E(£fa,tf). 

Similar to ( 13. 17ft . one needs to check whether the constants are non-zero: 

(7.18) (\fi\ Pl - 1 f 1,^)^0 and (£/,,#) ^0, 

which is not an easy problem and can lead to some restrictions for such a behaviour, 
though is crucial for any hope to see a bifurcation point. 

Under the conditions (17.181) . the parameter b(e) in (17.111) for p « p\ is given by 

(7.19) b (e) ^ {ll ( pi - p) ]Mn) ( s l I = ai ), e = Pl - p, 7* = ^f^) • 

The direction of developing in p of each p;-branch and whether the bifurcation is sub- or 
supercritical depend on the sign on the coefficient 7$. This can be checked numerically 
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only, but, in general, we expect that ji > 0, so that these nonlinear bifurcations are 
subcritical and the p\- branches exist for p < pi. 

Overall, the above formal analysis detects a number of key assumptions, which are 
necessary for such a nonlinear bifurcation to occur at the critical exponents pi given by 
(17. 8p . Recall again that, for n = 0, a more rigorous justification of the corresponding 
linearised bifurcation analysis is done in jl], where a countable number of p-branches was 
shown to be originated at bifurcation points (17. 9p . 

Overall, we claim that, for any n > 0, the ODE (17. 3p . with proper setting as y — > ±oo, 



(7.20) admits not more than countable number of p-branches of solutions, 



which are originated at the critical exponents (17. 8p (no rigorous proof is available still). 

7.3. Numerical experiments for k = 1. We briefly attempt to find numerical solutions 
of the equation (17. 3p . As usual, we look at the lower-order case, k = 1. Once again, in 
order to remove the nonlinearity in the highest (third)-order operator, the substitution 
Y = \f\ n f is used. This yields the semilinear third-order equation 

(7.21) f=\Y\-^hY: Y"< + \Y\~^Y'y + ^ \Y\~^Y - \Yf^Y = 0. 

Due to the complexity of the equation, which remains to be of the third order, there is 
not much hope of solving this problem using a shooting method, since, in fact, we do not 
know in detail the "nonlinear bundle" as y — > +oo. As usual and as in the semilinear case 
n = 0, there is a difficulty in finding the correct boundary points for y > 0, such that the 
correct oscillations are found. 

Hence, we return to the BVP setting, trying to "optimise" and "minimise" the oscilla- 
tory bundle for ?/> 1. However, even using this approach, which was rather effectively 
implemented in the simpler semilinear case n = in [H § 6], it is difficult to produce 
reliable numerics, due to the highly nonlinear /oscillatory nature of the problem, even for 
smaller values of p and n. In addition, we must be careful when plotting, as we must 
avoid approaching any nonlinear bifurcation points in p given in (17.81) . which actually are 
not known explicitly. 

As a first example, in Figure [TU we present "almost converging" VSS profile for p = 5 
and small n = 0.1. This example is of a particular importance: its tail for y ~ 35 is clearly 
not symmetric about {Y = 0} and is positive. According to 01 § 6], such a solution is not 
a VSS profile in R, and actually corresponds to some (obscure) BVP setting on a bounded 
interval, which is of no interest here. We should avoid such solutions in the future, even 
if these have been obtained with a perfect convergence up to the tolerances ~ 10~ 3 . 

The next Figure [TS] shows a better converging with a good tail VSS profile with the 
same p = 5 and a larger n = 0.6. We see that, for such larger n, the tail for y ^> 1 gets 
smaller and keeps being symmetric. 

For larger n > 1, the numerics gets less reliable, though we have got a number of 
"almost converging" results. For instance, in Figure [T6| this is done for p = 10 and n = 1. 
The tail is now larger for y ~ 12, and remains rather symmetric. 
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Y(y) VSS: p=5, n=0.1 ; "almost" converging tail (40000 points, max. res.=0.029) 
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Figure 14. A VSS similarity solution Y(y) of the ODE (17311) for p = 5 
and n = 0.1. 

Further increasing n requires also increasing p. In the next Figures [T7] and [181 we show 
VSS profiles for p = 12, n = 2 and p = 16, n = 3, respectively. In the former one, Figure 
[T71 we present two different profiles (solid and dotted lines), showing that non-converging 
tails do not affect the convergence in the dominant positive part of the profiles. 

Finally, we must admit that, since the correct "oscillatory bundle behaviour" as y — > 
+oo is still poorly understood, we cannot control it by choosing proper "minimal nonlinear 
components" (a proper symmetry is not enough; see [U § 6.3]). Therefore, we must 
confess that, even having good enough numerical convergence, we cannot guarantee that 
the above numerical examples get into the countable family indicated in f)7.20p . In other 
words, without a full use of correct "minimal oscillatory bundle" as y — > +oo (meaning 
non-posing any condition at the singular end-point y = +oo), the family of VSS profiles 
becomes continuous. Then, for any fixed value of p > n + 1, "solutions" f(y) represent 
a continuous line, instead of a predicted at most countable subset of points, lying on the 
above p-branches. 

On the other hand, it is plausible that those Figures correctly describe a general geom- 
etry of VSS profiles. In particular, we have always observed a strong "stability" of the 
first maximal positive hump, which turned out to be rather independent of the resulting 
tail for y ^> 1. In other words, we do note that, whilst the oscillatory part is difficult 
to obtain (and properly justify mathematically), the non-oscillatory structure is rather 
stable and does not change much, when changing the length in which the problem is 
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VSS: p=5, n=0.6; "almost" converging tail (100000 points, max. res. =0.037) 
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Figure 15. A VSS similarity solution Y(y) of the ODE (17311) for p = 5 
and n = 0.6. 

evaluated on. In particular, the numerical profile in the last Figure [18] was obtained with 
the worst maximal residual = 10.23 (actually meaning no convergence at all). However, 
we guarantee that this non-convergence takes place in the tail only, while the dominant 
positive part remains stable and "rather convergent", so we do not hesitate to present 
such a Figure here. However, in almost all other Figures, the convergence is much better 
and is not worse than 2 — 5 %. It may be also noted that this confirms such oscillatory (at 
least, symmetric) tails, even when these are large in the amplitude, are essentially zero in 
some "weak" sense (which we do not want to specify at the moment). 

7.4. A nonlinear limit n — > +oo: an example. Finally, let us note that, as in Section 
IU one way of finding proper oscillatory patterns f(y), would be to look at the behaviour 
as n — > +oo and reduce the ODE to a simpler one. We present an example of such a limit 
along the following straight line on the {n,p}-plane: 

(7.22) p = n + 4-»+oo as n -> +oo. 
Then the ODE (17311) reads 

(7.23) Y"> + ^ (\Y\-^kYy)> - \Y\&Y = 0, 
so that, on scaling, one gets 

n _i_ i ~ ~ ~ n ~ o~3~ 

(7.24) Y = (n + 3)- — Y Y'" + (\Y\~~Yy)' — (n + 3)~« = 0. 
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VSS: p=10, n=1 ; "almost" converging tail (60000 points, max. res. =0.054) 
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Figure 16. A VSS similarity solution Y(y) of the ODE fTf^TT) for p = 10 
and n = 1. 

Passing to the limit n — >■ +oo, in the ODE in (I7.24p . in the class of uniformly bounded 
solutions and using that (n + 3)~~ —> 1, yield two terms as in the equation (14 .4p with an 
extra linear one: 

(7.25) Y'" + (\Y\~^Yy)' -Y = 0. 

However, since unlike (I4.4p . this is a third-order ODE, an algebraic treatment of the first 
nonlinear eigenfunction Y (y), as in Theorem 14. 1[ becomes rather illusive. Anyway, this 
shows a principal possibility to study the "nonlinear" limits n — > +oo. 
On the other hand, scaling also the independent variable y, 

(7.26) Y = CY, y = ay, and C = (n + 3)~ 2 ^ i a E ^ li , 
yields, instead of (I7.24p . 

(7.27) Y'" + (\Y\-^Yy)' - (n + 3)"» a 3+ - \Y\^Y = 0. 

Therefore, passing to the limit n — > +oo in (I7.27p . after integrating once, we arrive 
precisely at the equation (14. 4p provided that 

a = a(n) — > + as n — >■ +cxo. 

Therefore, we obtain the same nonlinear eigenf unctions for / = 0,1,2 via the above 
algebraic-geometric approach, but, according to (I7.26p . on expanding subsets in the inde- 
pendent ^-variable. 
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VSS: p=12, n=2; example of non-converging tails and conv. positive part 
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Figure 17. A VSS similarity solution Y(y) of the ODE fl7T2TT) for p = 12 
and n = 2. 

Recall that the above limit as n — > +oo, with a possible study of branching as in 
Section 14.41 occurs along the straight line (17.221) (or in its "small neighbourhood" ) in the 
2D parametric {n,p}-plane. Along other lines, the limits can be different and lead to 
other patterns Yi(y), Y 2 (y), etc. 

Finally, overall, as we have seen, the ODE (I7.2ip represents a serious theoretical chal- 
lenge with respect to both analytical study (n-branching, p-bifurcation diagrams, and 
p-branches; the latter are known for n = [H § 7], etc.), as well as even a numerical one. 
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